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Motivated by the recent experiment on p-orbital band bosons in optical lattices, we study the- 
oretically the quantum phases of Mott insulator and superfluidity in two-dimensions. The system 
features a novel superfluid phase with transversely staggered orbital current at weak interaction, 
and a Mott insulator phase with antiferro-orbital order at strong coupling and commensurate flUing. 
We go beyond mean field theory and derive from a microscopic model an effective action that is 
capable of describing both the p-band Mott insulating and superfluid phases in strong coupling. We 
further calculate the excitation spectra near the quantum critical point and find two gapless modes 
away from the tip of the Mott lobe but four gapless modes at the tip. Our effective theory reveals 
how the phase coherence peak builds up in the Mott regime when approaching the critical point. 
We also discuss the finite temperature phase transition of p-band superfiuidity. 

PACS numbers: 03.75.-b, 05.30.Rt, 64.60.My, 67.85.-d 



I. INTRODUCTION 

Recent years have witnessed a dramatic experimental 
progress on cold atom systems in optical lattices P, Q ■ 
Detailed studies of strongly correlated quantum systems 
are possible thanks to the unprecedented controllabil- 
ity of cold atom experiments. The lowest band Mott- 
superfluid transition [3, 0] for lattice bosons has been 
successfully observed in experiments Experiments 
on populating bosons on excited bands are also put for- 
ward (y, [A- And a more recent p-band boson exper- 
iment [Sj, |9[ opens up a new thrust towards observing 
the exotic phases of bosons on higher bands with long 
life time. Given growing experimental progress on ex- 
cited band bosons, probing detailed features of p-band 
Mott insulators and p-band superfluidity is attracting 
broader interests. There are numerous theoretical work 
on excited band bosons focusing on proposing exotic 
phases [iMl, 

which have demonstrated the fascinat- 
ing physics associated with bosons on excited bands of 
optical lattices. The quantum phase transition from p- 
band Mott insulator to superfluid phase has been studied 
within the Gutzwiller mean field approach [l^ [l5| . For 
large interactions, the energy is minimized by an incom- 
pressible state with an orbital order. And for weak inter- 
actions the kinetic energy dominates over the interaction 
and drives the system into a superfluid with a feature of 
transversely staggered orbital current (TSOC) [ill- The 
competition between the kinetic energy and interaction 
energy is well described within the Gutzwiller approach, 
however the single particle spectra and the momentum 
distribution are out of reach within this approach. 

In this paper we apply the method of effective action 
beyond the Gutzwiller mean field, and explore the sin- 
gle particle spectra in both of the p-band Mott insulator 
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phase and the TSOC superfluid phase in two dimensions. 
We have studied the phase coherence in the Mott insula- 
tor phase and found that sharp peaks rise at finite mo- 
menta ((±7r, 0) for the p^; band and (0, ±7r) for the Py 
band) when the Mott gap is small. This offers new ap- 
proaches of preparing coherent matter waves from Mott 
insulators. From the p-band Mott insulator phase to the 
TSOC superfluid phase, the global U(l) symmetry and 
the time reversal T symmetry are broken. Away from the 
tip of the Mott lobe, we find two gapless modes at the 
critical point; while at the tip, we find four gapless modes 
due to the particle-hole symmetry. For TSOC superfluid 
phase we go beyond previous study in the weak coupling 
limit [ll|, [l^l and consider the leading effect of Hubbard 
interaction in the strong coupling regime. Our theory is 
capable of capturing the main feature of TSOC super- 
fluid phase in the strong coupling regime, where the crit- 
ical point of Mott-superfluid transition is located. The 
isotropy of the sound velocity of the TSOC superfluid 
phase is explained. Finally, the finite temperature phase 
transitions of the strongly interacting TSOC superfluid 
phase are discussed. 

The paper is organized as follows. In Sec. II, we dis- 
cuss the microscopic model describing the extended Bose- 
Hubbard model with p-orbital degrees of freedom and re- 
visit the features of the p-band Mott insulator phase and 
the TSOC superfluid phase. In Sec. Ill, we derive an 
effective action for the extended Bose-Hubbard model in 
the strong-coupling limit by performing two successive 
Hubbard-Stratonovich transformations of the inter-site 
hopping term. In Sec. IV, we calculate the single particle 
spectra and the momentum distribution for the p-band 
Mott insulator phase. We apply Bogoliubov theory on 
our effective action and study the sound velocity in the 
TSOC superfluid phase. The finite temperature phase 
transitions of TSOC superfluid phase are also discussed 
in Sec. IV. We conclude with a brief discussion in Sec. 
V. 
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II. MODEL AND PHASE DIAGRAM 

We start with a microscopic extended Bose-Hubbard 
model with p-orbital degrees of freedom on a square lat- 
tice [13, [m 

H = Ht + i?onsitc, 

Ht = J2r K(i")a2-(r + x) + al{r)ay{r + y) + h.c] 
-t± [al{r)ar,{r + y) + al{r)ay{r + x) + h.c] , (1) 
i?o„sito = Er ^iHr){n{r) - |) - ^L^ir^) ~ y.n. (2) 

Here, aj,(r) and a]^(r) are bosonic creation operators of 
Pa; and Pj^ orbitals at r. The discrete variable r labels the 
sites of a square lattice. The lattice constant a is set to be 
1. i (tj^) is the longitudinal (transverse) hopping between 
nearest neighbor sites, U the on-site repulsion and the 
local angular momentum operator = —ia^yUx + ial,ay. 
The average occupation number n of bosons per site is 
fixed by the chemical potential /i. 

Because < < and t± > 0, both of p^; band and p^ 
band show minima at finite momenta. It is inconvenient 
to take the long wavelength limit of the original lattice 
boson fields. To overcome this inconvenience we intro- 
duce the following staggered transformation 
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(3) 

(r) and V'J (r) are lattice field operators for pseudospin 
I t (r)> = i-r\p.)+i{-)y\Py) and | i (r)) = (-)-b.) - 
i{—)y\py) components, where \px) and \py) are local p^; 
and Py orbital states. 

In this pseudospin representation, the Hamiltonian 
reads 



H 



E 

r,r' 



with 



T(x) 



-n(r)--i,(r)^)-Ain(r), (4) 



2 

t-t 



(5) 



where a =t,l, "(r) = i'^{r)ijj^{r) + ?/'j^(r)?/'^(r), and 

L,{r) = (-l)^+y[i/;|(r)^/>^(r) - V|(r)V';(r)]. The fourier 
transform of T^^i gives 



e(k) 



en(k) en(k) 



(6) 



where e^^(k) = e^^(k) = {t - tj_){cos{kx) + cos{ky)) and 
en(k) = e;T(k) = (i + i±)(cos(fca.) - cos(fcj^)). 
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FIG. 1: (Color online). The phase diagram determined by 
the effective action method. The filling factor u of the Mott 
regime shown is 1. t± is set to be O.lt. TSOC means the 
transversely staggered orbital current superfluidity The 
alternating p^-Py pattern shown in the Mott regime is the 
pattern of the Mott insulator with filling of i/ = 1 in the p- 
bands. The staggered p^ ± ipy pattern in the TSOC regime 
illustrates the orbital current order in the TSOC phase. 



It can be verified that the band structure given by e(k) 
shows a minimum at zero momentum. Thus we have ob- 
tained a theory which is convenient for us to take the con- 
tinuum limit. The internal symmetry group of the Hamil- 
tonian is U(l) X T, where T denotes the time reversal 
symmetry (i/'f (;) (x, f) ^ i/'^(^)(x, -i)) and U(l) denotes 
the global phase rotation symmetry (tpa- e^^V'o-)- 

In the strong coupling limit, the system is in a Mott 
insulating phase for commensurate filling (filling factor v 
is an integer), where the filling is defined as the occupa- 
tion number of bosons loaded on p orbits per site. For 
filling factor v larger than 1, the interaction term favors 
local orbital current states because of the {—LD term 
in Hamiltonian (Eq. 2]); these vortex-like states form a 
vortex- antivortex pattern due to super-exchange [lol.[l5l|. 
For filling factor = 1, the local vortex- like states are no 
longer favorable because the interaction term does not 
contribute on single particle states. Mathematically, the 
operator term is equal to identity, when acting on the 
single particle states, and thus does not favor vortex-like 
states. The Mott phase has an antiferro-orbital order, 
i.e., an alternating p^-Py pattern (FIG.[T]), which breaks 
lattice translation symmetry [T^]. We are interested in 
the long wavelength modes within this phase. 

In the weak coupling limit, the system is in superfluid 
phase and the dispersion (obtained by diagonalizing e(k)) 
shows minima at zero momentum. The two minimal sin- 
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gle particle states carry lattice momentum k = and 
pseudospin a =t,i, and they are related by time rever- 
sal (T) transformation. Due to the (— i^) term in Hamil- 
tonian, the T symmetry is spontaneously broken in the 
ground state, i.e., either (-0^) or (-0^) is finite. It is clear 
from Eq. ([3]) that the original particles form a staggered 
Pa;±ipj, pattern (FIG. [T|) in this superfluid phase, which 
is named TSOC [ll|. Thus going from the Mott insu- 
lator phase with filling v ^ 1 to the TSOC superfluid 
phase, the U(l) x T symmetry is spontaneously broken. 
The phase transition (FI G. [U ) is confirmed by Gutzwiller 
mean field calculations [lO[. However, the momentum 
distribution and the correlation functions in the Mott 
phase are out of reach within Gutzwiller mean field cal- 
culation. Motivated by this, we develop a theory valid 
in the strong coupling regime. With this theory we cal- 
culate the single particle spectrum for both of the Mott 
phase and the TSOC superfiuid phase and discuss how 
the Mott gap closes at the critical point and how the 
phase coherence peak develops in the Mott phase, and we 
explain the isotropy of the sound velocity of the TSOC 
superfluid phase in the strong coupling regime. 



III. EFFECTIVE ACTION 

To capture the main feature of p-band Mott insu- 
lator and the TSOC superfluid in the strong coupling 
regime, we aim at a theory capable of incorporating 
the local Mott gap, which is the leading effect of the 
Hubbard interaction, in a non-perturbative manner. To 
do this, we follow the procedure in Ref. [iGl - flsj . We 
first write the partition function Z as a functional in- 
tegral over complex fields ij^a with the action S[ip* ~ 

lo dTY.a,r {M-^)drM^) + Wc iutroducc au 

auxiliary field 0cr to decouple the inter-site hopping term 
by means of a Hubbard-Stratonovich transformation and 
obtain 

= Zo y" i?[C0,]e'^^^.".'"^"'(exp[(0|V)+c.c.])o 
= Zo y i?[C,'/'.]exp(CT-;,0,, +iy[0;,0,]), (7) 

where the shorthand notation {4'\''P) = 
Ecr r dT(j)*^{v)i}}c{r) and T^^ denotes the inverse of the 
hopping matrix. 5*0 and Zq are the action and partition 
function in the local limit (i,ij_ = 0). (. . .)o means av- 
eraging over the local action Sq , • The introduced 
generating function Vl^[0*,0cr] = ln(exp ((0|i/)) -I- c.c.))o. 
The local action Sq, which is equivalent to the original 
action without tunneling term, is invariant under a 



expansion of W^[(/)* , <f>„] respecting this symmetry yields 

W ^ J dTidT2^Ga{r,Ti - T2)(j)*„{r,Ti)(j)a{r,T2) 



Xa,a, (r, 1234)0;^ (1)0,, (2)0;^ (3)0,, (4) 
+O(06), (8) 

where the indices 1, 2, 3, and 4 indicate time ti, • • • , T4. 
And because of T (time reversal) symmetry, G^ — G^, 
and xtl- = X^^ • In the transformed theory, the quadratic 
term of in Eq. ([7|) is dominant in the strong coupling 
limit, t/U — > 0. We thus truncate the power expansion 
to quartic order. Now that the fluctuations of fields 
are controlled by (Eq. ([7])), perturbative renormal- 
ization group (RG) analysis finds higher order terms are 
irrelevant. From the definition, the coefficients are read- 
ily obtained as follows 

G,(r,ri -T2) = (■0,(r,ri)0';(r,T2))o, 
X.,,,(r;1234) = (7^,, (1)0^, (2)^^., (3)^. (4))C- (9) 

The effective action for fields is S'cff [0* , 0a] = 
/ dr Er,r' -C(r)?^<r; (r - r')0a' (r') - M/[0;, 0.], which 
can be used as a starting point to study the instability 
of Mott phase with respect to superflui dity by treating 
fields as superfluid order parameters [l^. However it 
is inconvenient to calculate the excitation spectrum and 
the momentum distribution from this action [l6| . Also, 
the theory S'off[0* , 0cr] does not provide a clear picture of 
quasi particles. The above difficulties can be overcome 
by performing a second Hubbard-Stratonovich transform 
following the method of Rcf. [l^ , 

f i?[<^>0*0]e-^^^-'^-'-[(^l«+=-^-l+'^K,^-l. 



(10) 

Integrating out fields gives the effective action 

s,s[^:,^,] = -w[^;,^,] + J dT £0,(11) 



C/(l) X U{1) transformation ipa- 



J0C 



^a- By definition. 



W\ 



ha-] is independent of hopping (t and tj_). Power 



<ri,r2> 

with 

-Co = tpl^{ri,T)Ta,a2iri - r2)v5a2(r2,r), 

ly = ln(exp(-[(^|0) + c.c.]))H', (12) 

where ( - f D('i't:'P^)i-)-MW[K..<p^]) 

wneie [...)w - f D{r^,4>.e^p{w[4>:,,4>.]) ' ^^^^^ 
the functional VF[0*,0ct] is independent of hopping, 
VF[(/3* , ifa] is also independent of hopping by definition. 

Now the task is to calculate the functional W. The 
essence is to evaluate the expectation value of the expo- 
nential operator in a system described by the effective 
action W. Note that the action W (Eq. (|5])) is a power 
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expansion of the small parameter l/U (strong coupling), 
so the perturbative field theoretical method is valid and 
powerful to compute the expectation value. In this man- 
ner, we perform the power expansion of W[ip'^,ipa-] re- 
specting the U{1) X C/(l) symmetry and obtain the effec- 
tive theory for ipa as follows 

/ (irifir2^(y9;(r,ri)G^^(ri - T2)^<,(r, rs) 

+ / X! ¥'ai(l"l,T)TlTi<T2(l"l - r2)V5<T2(l"2,T) 

<ri,r2> 

+ l9.ia,J dT^|(^,,(r,T)|2|<^,,(r,T)|2. (13) 

Here we have truncated the power expansion to quartic 
order. The vertex term can have a finite range behav- 
ior in T space in principle, but we can take the static 
limit because the r dependent corrections would be ir- 
relevant in the RG sense. The field i/j describes the 
bare (original) bosons with onsite interaction J7, while 
the field ip describes quasi particle excitations which are 
greatly suppressed by the energy gap (of the order oi U) 
in the strong coupling regime (this physical understand- 
ing comes from Eq. (jASp V Thus the quasi particles are 
dilute, and we expect the three body scattering process 
is negligible, which further guarantees the validity of the 
truncation performed in Eq. (|13p . Here, we want to 
emphasize the quasi particles hop around through T^^' 
(Eq. (|13p V The quasi particles are mobile instead of lo- 
calized. The dispersion of the quasi particles is discussed 
in the next section. 

Despite the difference of bare correlators of ip and ip, 
it is proved that the correlators defined by the original 
action S[tjj*,ip] are exactly equal to that defined by the 
infinite series of power expansion of this effective action 
Scff[ip*, ipa] (Appendix B). With valid truncations, the 
connected correlators of ip reproduce the connected cor- 
relators of ?/) approximately. In principle, one can follow 
the above procedure and get G"^^ (ri — T2) and 5ctict2 • K is 
straightforward to calculate these coefficients for the one 
component Bose-Hubbard Model. However it is inconve- 
nient to proceed in this approach due to the complexity 
induced by the local degeneracy of ground states of the 
action So- Since these coefficients are independent of hop- 
ping, we decide to calculate the coefficients by identifying 
correlators of this effective action and the original corre- 
lators (defined by the original action) in the local limit 
(t,t± = 0). The derivation and the results are summa- 
rized in Appendix A. 

After the above manipulation, the leading effect of 
the Hubbard interaction, namely, generating a local en- 
ergy gap in the quasi-particle spectrum, is included in 
the quadratic part of Scsi^p'^, ifia-], which makes the fol- 
lowing Bogoliubov analysis valid in the strong coupling 
regime. With the effective theory Scs[ip'^,ipa], the su- 
perfluid phase in the strong coupling regime is described 



as a superfluid phase of weakly interacting quasi parti- 
cles. We want to emphasize that this physical picture is 
lacking in the theory S'efi [</>*, (pa]- 



IV. ZERO AND FINITE TEMPERATURE 
PHASE TRANSITIONS 

In this paper, we focus on the lowest Mott lobe regime 
{1/ — 1) for which time reversal T symmetry is not broken. 
We start our analysis from the effective action (Eq. P^ ). 
We use the correlators of quasi particle fields ip, which 
can be calculated within Bogoliubov theory, to approx- 
imate the correlation functions and momentum distri- 
bution of the bare particles. (The reason is explained 
in Appendix B.) The method, we apply here, yields re- 
sults jlal that agree with for example the RPA calcula- 
tion |18j . when calculating the momentum distribution 
of the s-band Bose-Hubbard model. 

First, let us determine the ground state which is ho- 
mogeneous and static after the staggered transformation 
(Eq. ([21)). The free energy functional is 

= (^CT Ho) + e<T<T(k = Q))ns,a + 55^1^2 "s,<ti «s,a2 , 

. (^^^ 

where is the number of lattice sites, Qa{iijj) is the 
fouricr transform of G(t(t) and Us.a = \{^a)\'^ is the 
superfluid density of the pseudospin a component. By 
minimizing the free energy functional, we get 

^^^^(Q^lm^^i^ if^.-HO) + ...(k = 0)<0; 
1 otherwise. 

(15) 

In the Mott regime, n^ cr vanishes. In the superfluid 
regime, only one component is finite, i.e., either Ug^^ 
or n^^^ is finite and the other vanishes, because the off- 
diagonal part g^^ is greater than the diagonal part g^^ 
{= g^i). Thus U{\) X T is spontaneously broken across 
the Mott-superfluid phase transition in this model for fill- 
ing V = 1. The phase boundary is shown in FIG. [1] The 
phase diagram is consistent with Ref. The Mott 

regime determined by our approach is larger. (When 
applied to calculate the phase boundary of s-band Bose- 
Hubbard model, the method adopted here yiel ds results 
that agree with other mean field theories [la)-) In the 
following part, we assume the superfluid component is 
the 't' component. The superfluid density scales as 
ns^-\ ^ I?! — away from the Mott tip when the mean 
particle number per site, n, is not equal to v, whereas it 
scales as n^,^ ~ (t — tc)^ at the Mott tip where n = v- 

Next, we explore the fluctuations (pai'^i t) ~ 'y3cr(r, t) — 
yJUs.rj- Expanding the action to the quadratic order of 
the fluctuation flelds 'pa, we get 

S\(pl,q^a\ = i^vi/t(k,ic.)[/i(k,zw)]M'(k,zL.), (16) 

k,Lj 
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with 



(17) 



2 



and 



hii 
his 

/l44 



^21 = ffrr^s , 
^31 = en(i<^) ' 

+ ^^^^ (tw) + .g^;n^ , 



(18) 



where (^(T(k, iuj) is the fourier transform of <^cr(r, t). Only 
non-zero [h] matrix elements are listed above. Using the 
quadratic action for fluctuations, we calculate the Bogoli- 
ubov spectrum. 



A. Mott phase 

In this part, we study the momentum distribution of 
p-band Mott insulator phase and we show that diverg- 
ing peaks for bosons on p^, and Py bands rise at finite 
momenta when the Mott gap closes. This offers new 
approaches of preparing coherent matter waves in exper- 
iments. 

In the Mott phase, superfluid density Ug = 0. The 
Green function ^#ctict2 (k, it^) = (V'cti (k, iw)'!/'CT2 (^i ^i^)) 
(k, iw)<p*2 (k, iw)) is readily obtained as 
e^jCT2(k) -I- g^l{iuj)5„^a^. Solving 
the equation det[^#^^(k, w)] = 0, we get the single 
particle spectra shown in FIG. [5] and FIG. [31 Deep 
in the Mott regime, all single particle excitations are 
fully gapped. It can be verified that the Green function 
?#CTi,cr2 (k, «i^) is diagonal in the {■0a;(k, iw), ij}y{k,iu!y\ 
basis. Thus we can label the spectra by and . 
The E^{E^) branch is the p^, band of particle (hole) 
excitations; while E^{E~) branch is the Py band of 
particle (hole) excitations. Upon the phase transition 
point, the gap at k = drops and approaches zero. 
Away from the Mott tip, two particle branches (i?+ 
and £■+) touch zero first (FIG. [3]), which causes the 
instability of the Mott insulator phase and drives the 
phase transition from Mott insulator to TSOC super- 
fluid. For a more physical case in which the density 
is flxed to be an integer, both of particle and hole 
branches close simultaneously when increasing hopping 
because of the particle-hole symmetry at the Mott tip 
(FIG. [3]). For both cases the gap closes right at the 
phase transition point determined by minimizing the 
free energy functional in Eq. (jlSp . 

At zero temperature, we obtain the momentum dis- 
tribution n(k) EE Tr(?/Ci (k)^/v2 (k)) of the p-band Mott 
insulator phase from the Green function as n{k) — 
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FIG. 2: (Color online), (a) The single particle spectra along 
the fcj; axis (shown in the inset) deep in the Mott regime 
with parameters fj./U = OA,t/U = 0.02, fx = O.lt. (b) The 
single particle spectra along the axis (shown in the inset) 
in the Mott regime near the critical point with parameters 
/i/C/ = 0.4, t/U = 0.065, t± = O.lt. The Mott gap drops when 
increasing hopping and the gap of particle branches closes at 
the critical point for fi/U — 0.4. 
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FIG. 3: (a) The single particle spectra along kx axis deep 
in the Mott regime with parameters ^/U = 0.27, t/U — 
0.02, fx = O.lt. (b) The single particle spectra along the k^ 
axis in the Mott regime near the critical point with param- 
eters fj./U = 0.27, t/U = 0.09, tx = O.li. The Mott gap for 
both of the particle and hole branches closes at the critical 
point at the Mott tip regime. 



^ J c?a;e*"*'^^cria2 (k, iw). The momentum distribution 
measures the spectral weight of the negative pole of 
Green function [^^(k, w)] (the negative pole refers to the 
real part of the energy spectrum being negative). Deep in 
the Mott regime, the momentum distribution is very flat 
(FIG.|4(ay I. Approaching the phase transition point from 
Mott insulator side, a peak at zero momentum develops, 
and the peak diverges right at the critical point where the 
Mott gap closes(FIG. |4(d)| . The continuous development 
of the peak of the momentum distribution at k = means 
the phase coherence {{ipl.{r)'ipx{'r')) and (?/'J(r)-0.y(r'))) 
develops continuously when increasing hopping in the 
Mott regime. However, we do not see the phase coherence 
of ipx and ipy components, i.e., {ipl.{r)ipy{r')) always van- 
ishes in the Mott regime even when r = r'. To compare 
with experiments, we also obtain the Green function for 
the original bosons and get the momentum distribution 
of original bosons shown in FIG. [5] Near phase transition 
point in the p-band Mott insulator phase, the momentum 
distribution for p^^ bosons shows a peak at (±7r, 0) and 
the momentum distribution for py bosons shows a peak 
at (0, ±7r). This offers possibilities of observing coherent 
matter waves from p-band Mott insulator. 
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FIG. 4: (Color online). The momentum distribution defined 
by Tr(-i/)Jj (k)'!/'a-2 (k)) in the Mott regime with parameters 
^j,/U = 0.3, tx = O.lt. For (a) through (d), the parallel hop- 
ping t/U is 0.02, 0.04, 0.06 and 0.09 respectively. A coherent 
peak rises continuously when approaching the Mott-superfluid 
transition point from the Mott insulator side. The unit for 
lattice momentum k is with a the lattice constant. 
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FIG. 5: (Color online). The momentum distribution (a) for 
Px band boson and (b) for py band boson in the p-band 
Mott insulator phase near the Mott-superfluid phase transi- 
tion point. The parameters we use are fi/U — 0.3, t/U = 0.09 
and t± = O.lt. The unit for lattice momentum k is a^^ with 
a the lattice constant. 



B. Superfluid phase 

The single particle excitation spectrum of the TSOC 
superfluid is calculated within Bogoliubov theory. Solv- 
ing the equation det[h{k, iuj)] = 0, the spectrum is ob- 
tained. Deep in the TSOC superfluid phase, there are 
four branches, one of which is gapless and linear around 
momentum k = (FIG. |6]). The spectra can no longer 
be understood in the same way as in the Mott regime, 
since the time reversal symmetry is broken. For momen- 
tum kx = iky, the excitations have definite pseudospin, 
while for generic momenta, the excitations do not carry 
definite pseudospin. The ^ and I components are actu- 
ally mixed unless k^ = ±ky. The physical reason for this 
is that the dispersion 6^1^20^) is not diagonal for generic 
momenta. 

Deep in the superfluid regime, the off-diagonal term 
e^l does not affect the U(l) phase mode for the reason 
that the mode of flipping pseudospin is fully gapped (the 
gap is g^iUg) in that regime (FIG.[6jb)). And the spec- 
trum of the phase mode (i^;7(i) (k)) of TSOC superfluid is 
fully determined by hu, hi2, /121 and /i22 in Eq. ([T8| . all 
of which are isotropic to the quadratic order of momcn- 



>.2 


2 


(a) 





1 2 3 

k (units of 1/a) 




12 3 

k (unitH of 1/a) 



FIG. 6: (Color online). The Bogoliubov spectra in the TSOC 
superfluid phase, (a) shows the spectra along k^ axis near 
critical point with parameters n/U — 0.5, t/U — 0.065, t± = 
O.lt. (b) shows the spectra along kx axis deep in the superfluid 
phase with parameters fi/U — 0.5, t/U = 0.1, t_L ~ O.lt. 



tum k. The sound velocity defined by 9kWc/(i)(k)|k-i.o is 
thus isotropic in the weak coupling regime as derived in 
Ref. In the strong coupling regime, the gap of flip- 
ping pseudospin becomes smaller near the phase transi- 
tion point (FIG.iniJa)), and the coupling (e-|-^(k)) between 
-tp-f and ipi 'v^ill modify (^[/(ij (k) . However, the correc- 
tion to W(7(i)(k) in the limit of k — > is to the order of 



2iK-K) ^ 



(k)|i 



j.0 is unaffected. Thus 
the sound velocity of TSOC superfluid is isotropic in the 
strong coupling regime. 



C. Finite temperature phase transitions of TSOC 
superfluid phase. 

In the TSOC superfluid phase, we can neglect the tem- 
poral fluctuations of the fields. Thus we have 



vW2 



J9(t) + i 



i-y 
(-) 



(19) 



where s^/y(r) = or 1. This substitution captures the 
thermal fluctuations of the TSOC superfluid phase in the 
strong coupling regime. And the energy functional de- 
scribing the fiuctuations of 6 and s fields is 

E = -Jg [cos(e'(r) - e{r + x)) + cos(6'(r) - e{r + y))] 

r 

|^_ J^g*7rsy(r)g?7rsy(r+i/) _ j^^i7ZSy{r) ^J7rsy{r^x) 

r 
r 

where Jg = {\t\ + \t±\)ns/2, Ji = \t\ns/2 and J2 = 
\t±\ns/2. The 9 part is an isotropic XY model and 
the Sx Sy parts arc anisotropic Ising model. Consid- 
ering the thermal fiuctuations of 9{r), Sx{r) and Sy{r) 
separately, we can define two transition temperatures — 
the Kostcrlitz Thoulcss (KT) transition temperature 
Tkt and the Ising transition temperature Tjsing. Tkt 



is well approximated by Tkt 
Ising temperature is Tiding 



The 



-Jeff, where JeS 



log(l+V2)' 

is determined by sinh(2 Ji/Tising) sinh(2 J2/Tising) 
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\ > 



T=0 
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Appendix A: The coefficients of tlie effective action 



FIG. 7: (Color online). The schematic plot of finite temper- 
ature phase diagram in the TSOC superfluid regime. Tci is 
the Ising transition temperature. Below Tci, the orbital or- 
der and the superfluid order coexist for the TSOC superfluid 
phase. Tc2 is the KT transition temperature. For temper- 
ature above Tci and below Tc2, the orbital order no longer 
exists but the superfluid order still survives. Above Tc2, no 
order exists and the lattice bose gas is in the normal phase. 



sinh^ (2 Jeff /Tiding) |2ll. For \t^\ < Ti,i„g « Tkt be- 
cause Jeff — ^ when tj_ — > 0. At temperature lower than 
Tising, the TSOC superfluid phase has an algebraic U(l) 
ordering and a long range Ising ordering (orbital order). 
At temperature higher than Tiging, the orbital ordering 
disappears. Because of the existence of free Ising kinks 
(flips of pseudospin), half vortices with boundary con- 
dition § dO = i: are possible, which will modify the KT 
transition temperature. The modified KT transition tem- 
perature is Tkt ~ '^Je when <C The schematic 
finite temperature phase diagram is shown in FIG. [T] 



V. CONCLUSION 

In this paper we have developed a theory of an ex- 
tended Bose-Hubbard model with p-orbital degrees of 
freedom, which is valid in the strong coupling regime. 
With this theory, we explored the single particle spectra 
in the Mott insulator phase and the Bogoliubov quasi- 
particlc spectra in the TSOC superfluid phase. We stud- 
ied how the momentum distribution develops in the Mott 
insulator phase when increasing the inter-site hopping 
and found that diverging peaks rise in momentum dis- 
tribution at finite momenta when the Mott gap closes. 
We explained the isotropy of the sound velocity of the 
TSOC superfluid phase. The finite temperature phase 
transitions of the strongly interacting TSOC superfluid 
phase are discussed. 
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We calculate the coefficients of the effective action 
in Eq. ([T5|) by identifying the connected correlators of 
the effective action and that of the original Hamiltonian 
(Eq. |4]). The connected correlators of the original local 
action is calculated in the operator representation in the 
occupation number basis. 



1 



Gt(r,T) = — Tr[T,^t(i-,T)^j(r,0)e-'^^°] 



^^(7iTO|e-('3^")^>t(r)e-"^''V':f(r,0)|nm), 

]^^g-(/3--r)«(n.ni)g-Te(n-|-l,m) (Al) 



-E 



where \nm) = 



-|0) are the eigenbasis of 



the local interaction Hq with eigenvalues e(n, m) = 
^ ((n -I- m)^ - |(n + m) - i(n - m)^) - ^{n + m). The 
fouricr transform of this correlator gives 

g^itio) = f o^TG^(r)e'"^ 



1 f g-/3e(n+l,m) 
— ^(n + 1) —— ^ ^— ^ 

Zq ^-^ I tuj + e(n, m) — e(n + 1, m) 

-/3£(n,in) 



iuj + e(n, m) — e{n + 1, m) 



(A2) 



In the low temperature limit (/?[/ 3> 1), the exponen- 
tial term selects out the local ground state, contributions 
from other states being suppressed. The double degener- 
ate ground states are |1) = \nQ,0) and |2) = |0, uq), where 
uq is defined by minimizing e{n,m). These two ground 
states are related by time reversal symmetry. Because of 
the double degeneracy, we further define correlators with 
respect to one single ground state. Q^'^ {iuj) is defined 

corresponding to the ground state |1), while C/| is 
defined corresponding to the ground state |2). Thus, 



no 



iuj + e{no — 1, 0) — e(no, 0) 
no + 1 



iuj + e(no, 0) — e{no + 1,0)' 
gi '(IUJ) = - — 

I 7/ , 



iuj + e{0,no) - e(l, ^q) ' 



(A3) 



It can be verified that g^^"^ = g'^^ and g'^^^ ~ ^j^'- Up 
to this point, the quadratic part of the effective action 
in Eq. is obtained. In the following, we always split 
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the correlators into two parts (labeled by superindices 
1,2) according to two ground states |1) and |2). In this 
paper we focus on the Mott insulator with filling ly = 1 
for which the time reversal symmetry is not broken. For 
the vortex-antivortex Mott insulator with filling v > 1, 
the local ground state will spontaneously choose either 
|1) or |2). And thus one can calculate all the correlators 
assuming that the ground state is |1) or |2) instead of 
taking the average. 

The four point function is given as 



Similarly the off-diagonal part is calculated as follows 
1 



-(1,2) 



Tri,2 [e~^"'>Trie/'"°ij^e~^'"°) 



~T3i/o\„/,t 



(A9) 



where Tri^2 means taking the trace with respect to the 
ground state |1) or |2). After some straightforward cal- 
culation, 



^0 



(AlO) 



(A4) 



To calculate the coefficients (7o-i<t2 in the static limit, we 
are only interested in the time average of x- 



Xo 



J dTldT2dT3X<yi<T2{Tl,T2,T3,0) 



Xflf2 ~l~ XlTi(T2 



(A5) 



The diagonal part is calculated as 
1 



-(1>2) 



Tri.2 [e-'^"°Trie-''"°'iP-^e-^'"'>) 



(1/2), 



(A6) 



where Tri.2 means taking the trace with respect to the 
ground state |1) or |2). After somewhat tedious calcula- 
tion we get 



-4(no + l)(no + 2) 



;(no, 0) - e{no + 1, 0)]2[e(no, 0) - e(no + 2, 0)] 
-4(no - l)no 



eK, 0) - e{n - 1, 0)]2[e(no, 0) - e(no - 2, 0)] 
-4no(no -|- 1) 



e(no - 1, 0) - e(no, 0)]2[e(no + 1,0)- e(no, 0)] 
4-4no(no + 1) 



;(no, 0) - e(no + 1, 0)]2[e(no, 0) - e(no - 1, 0)] 



;(no-l,0)-e(no,0)]3 
4(no + 1)2 



-(2) 



e(no,0)-e(no + l,0)]3 

-8 



(A7) 



[e(0,no) - e(l,no)]2[e(0,no) - e(2,no)] 
4 



[e(0,no)-e(l,no)]3 



(A8) 



Because of time reversal symmetry other four point cor- 
relators are readily obtained by x|^| = xj^j? xfi = X^\- 



with 
fo - 



-1 



and 
/i 



;(7io, 0) - e{no - 1, 0)]2[e(7io, 0) - e(no - 1, 1)] 
1 



[e(7io,0) 


-e(7io-l,0)]2[e(no,0)-e(no,l)] 
1 


[e("o,0) 


-e(7^o,l)]2[e(no,0)-e(?^o-l,0)] 
-1 


[e(no,0) 


-e(no,l)]2[e(no,0)-e(no-l,l)] 
-2/[e(no,0)-e(no-l,l)] 


[e(rto, 0) 


-e(7io-l,0)][e(no,0)-e(no,l)]' 
1 


[e(no,0) 


-e(no,l)]2[e(no,0)-e(no + l,0)] 
1 


[e(no,0) 


-e(no-f l,0)]2[e(no,0)-e(no,l)] 
-1 


[e(«o,0) 


-e(no,l)]2[e(no,0)-e(no + l,l)] 
-1 


[e("o,0) 


- eino + 1, 0)]2[e(no, 0) - e(no + 1 
-2/[e(no,0)-e(no + l,l)] 



(All) 



■{no, 0) - e(?io + 1, O)][e(no, 0) - e(7io, 1)] 



(A12) 



is obtained from 



with e{n,m) substituted by 



-ii) 



and 



-(2) 

t{m,n). For rip > 1, one can safely let x^^ — x|^^ 

thus = X^j,'' because e(n,m) = e(rn,n). However for 

no = I, both x^^-* and x"^^ are singular because e(l, 0) = 
e(0, 1); while the average of these two is finite by taking 
the proper limit e(l, 0) = e(0, 1) -I- 0+. 

Up to this point, we have obtained the four point cor- 
relators from the original Hamiltonian in Eq. In 
order to calculate gcna^-, we still need the time average of 
the four point correlator defined by the effective action in 
Eq. (fT3|) . Since the local ground state of this theory is not 
unique, one should not naively apply the Feynman rules 
of the usual 0^ field theory. Similar to the approach we 
used above, we calculate the correlators on each ground 



9 



state and then take the average. Thus we have 

{t1,T2,T3,0) 



Al) , ^(2) 



+ [(^(f (0)C?(2'(0))' + {Q?} f5„^...]} /2. (A13) 
Then gaia2 is obtained, 

+ [(^^^'(0)t^(',H0))2 + (t;(^))4j^^,,j}. (A14) 
We verify that these coefBeicnts satisfy g-^^ > g^^. 



Appendix B: Validity of the double 
Hubbard-Stratonovich transformation 



The connected correlators of the original boson fields 
i/'cr can be obtained from the generating functional 

Z[J*,J] = e~i'aTaa'iJa' " So[tp*,tp] + [( J|l/') + C.C.lj-g^) 

For example, 

where x = (r,T) and (. . ■)s[0*,i/'] means an average de- 
fined by the action Slip* = ip*Tca'il^a' + So ['/'*, V'] ■ 
Introducing the first Hubbard-Stratonovich transforma- 
tion, 

z[j*,j] = J D[r,i^; 0*, 0] exp {r.T^u-' - Som.i^a] 

+ [(0|^)+c.c.] + [(J|V') + c.c.]}. (B3) 



After a shift (per ^ 'Pa ~ Ja, we get 

Z[J*,J] =/£>[?/'*, V','/'*,0]e[('^l'^)+^-=-l-'^°['^J''^"l 

X exp { (0, - JaYT-;, {cP,. ~Ja')]. (B4) 

After integrating out the -0^ fields, we get 

Z[J\J] ^ZoJ D[(P*,cP] exp {(0, - j„)*t;^,{4>„, - J,0 

(B5) 

Now we introduce the second Hubbard-Stratonovich 
transformation 

Z[J\J] 

= ^0 ! D[ip* ,ip;(j)* ,ct)]e-xp{-f*aTc,c,''^a' + [{ip\J) + c.c]} 
X exp {-[(^10) + c.c] + W[cP%,ct,,\} . (B6) 

Integrating out the (per fields, we get 

Z[J*,J] = Zq I D\i^* ,ip\G^^^-(p*jr„„'ipa' + ^VP*a,^a\ 
-f[(^|J)+C.C.]}. ^ (B7) 

The generating functional for (^o- fields is equal to the 
generating functional of ip^ (up to a constant indepen- 
dent of the sources J). And this proves the connected 
correlators of ip„ are the same as that of ipa, although 
the action of (p^ is very different from that of ipn. For 
example, 



(:ipa{x)i}l,{x'))s{^*.^^\. 



(B8) 



where (...) s^hYp* ,ip\ means an average defined by the ac- 
tion S'eff[<^*,</?] = (^;T„a'<y3CT' - W{ipl,<^a\- 
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